Crystallization dynamics of a single layer complex plasma 
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We report a series of complex (dusty) plasma experiments, aimed at the study of the detailed 
time evolution of the re-crystallisation process following a rapid quench of a two dimensional dust 
liquid. The experiments were accompanied by large-scale (million particle) molecular dynamics 
simulations, assuming Yukawa type inter-particle interaction. Both experiment and simulation 
show a. Oct" (power law) dependence of the linear crystallite domain size as measured by the bond- 
order correlation length, translational correlation length, dislocation (defect) density, and a direct 
size measurement algorithm. The results show two stages of order formation: on short time-scales 
individual particle motion dominates; this is a fast process characterized by a = 0.93 ± 0.1. At 
longer time-scales, small crystallites undergo collective rearrangement, merging into bigger ones, 
resulting in a smaller exponent ct = 0.38 ± 0.06. 



From the very beginning of laboratory complex (dusty) 
plasma research, one of the main motivating and promis- 
ing features of these strongly coupled many-particle sys- 
tems has been the possibility for modeling classical collec- 
tive phenomena occurring in atomic matter on a size and 
time scale that allows direct observation at the particle 
level Experiments addressing e.g. transport phenom- 
ena in two dimensions (heat conductivity 0, viscosity 
Q, and self-diffusion H) dislocation dynamics 01, vari- 
ous melting scenarios [6|-|9| , and some aspects of freezing 
[lol [Tlj have already been carried out. 

The fact that a two-dimensional, hexagonal crystalline 
structure could easily develop in a laboratory complex 
plasma experiment was clear from the early days of the 
field [12, even though the theoretical background ex- 
plaining crystallization in low dimensions |13| is still in- 
complete and the subject of ongoing debate |14{ . 

Here we study the phenomena of pattern formation oc- 
curring when a liquid is rapidly quenched to a solid. Our 
qualitative expectation is that the initial amorphous liq- 
uid structure should evolve toward a ground state crystal 
through the process of domain coarsening, i.e., merging 
initially formed small crystallites into fewer and bigger 
ones. To characterize this process we measured the aver- 
age linear size of the crystallites and recorded their time 
evolution. Size measurement is based on: (i) the bond- 
angular correlation length, (ii) the translational correla- 
tion length, (iii) the inverse defect fraction, and (iv) a 
direct size measurement algorithm similar to the stan- 
dard "flood fiH" method used in simple graphical tools. 

Ref. [31 reported a time dependence of the orienta- 
tional correlation length in the form with a « 1/4 
measured in a two dimensional temperature quenched 
experiment in a single layer of spherical block copolymer 
microdomains in a thin film. Recent experiments with 
superparamagnetic particles (cx Ijr^ interaction) have 



shed new light on the dynamics of crystallite formation in 
colloidal systems resulting in a currently accepted value 
of a ~ 0.3 16j. Recent equilibrium melting studies, on 
the other hand, expose an important difference between 
colloidal suspensions and dusty plasmas: dipole systems 
(with cx l/r^ pair-potential) seem to follow the "dislo- 
cation unbinding melting" (KTHNY) picture jl^l, while 
in complex plasmas (with cx exp(— t/Ad)/?' interaction) 
rapid domain wall formation is observed, as described by 
the "grain-boundary melting" scenario [9| . 

Various numerical methods (including Monte Carlo, 



molecular dynamics, dislocation dynamics 17| and 
phase-field crystal simulations) were applied showing 
power-law typ e time evolution of the domain size, as re- 
ported in |l8| and (l9| with a « 1/3 for spinodal decom- 
position, in [20| with a ~ 1/3 for a polymer solution, 
in i21|] with a « 0.47 for a binary solution, in with 
a « 1/4 for block copolymers, and in 23] with a « 0.42 
for a nematic liquid crystal. 

We have performed a series of complex plasma ex- 
periments using the CASPER dust y p lasma experimen- 
tal setup introduced in detail in j24| . Mono-disperse 
melamine-formaldehyde micro-spheres with a diameter 
d = 6.5 /im ±1% were used. An argon gas discharge 
was driven using a P = 1.5 W radio-frequency source 
at 13.56 MHz at a pressure p = 4.7 Pa and a constant 
gas flow of about 10 seem. A single-layer of dust par- 
ticles was horizontally illuminated by a red laser diode 
and the scattered light was captured by a "FASTCAM- 
1024PCI model lOOK" fast CCD camera from the top 
at 125 frames per second with a resolution of 512 x 512 
pixels and a field of view (FOV) having a side length 
i = 30 mm. The vertical extension of the layer was mon- 
itored using a side-view camera. Our dust clouds con- 
sisted of approximately 2000 particles within the FOV. 
Without external excitation, the system rapidly settled 
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into a ground state configuration exhibiting three or four 
large crystal domains (forced by the confinement geom- 
etry of circular symmetry). The clouds were very sta- 
ble with no rotation or any other dynamical frustration 
observed. Melting of the single-layer was induced by 
temporarily introducing a ^ 20 Hz sinusoidal modula- 
tion to the DC self-bias of the powered lower electrode, 
causing a vertical shaking of the layer, which coupled 
to the horizontal motion of the particles (as studied in 
detail in [i^), resulting in the complete elimination of 
the crystalline order within the system. Rapidly switch- 
ing off the modulation allowed the system to settle into 
its vertical equilibrium position within less then '2,1'^^^^., 
where I'dust = 8.8 ± 0.4 Hz is the observed dominant 
dust particle oscillation frequency. Both the melting and 
crystallization cycles were repeated several times in the 
experiment and images of the particle suspensions were 
collected continuously. 

Figure [T] illustrates the principle data evaluation steps 
for two snapshots, one representative of the system in a 
liquid phase, the other in a solid configuration. First, the 
raw images were processed using a particle detection al- 
gorithm based on the center-of-mass method introduced 
and optimized in ^2^, providing accurate sub-pixel res- 
olution (Fig. [T](a)). Based on these position data, De- 
launay triangulation was then performed to find the true 
nearest neighbors (NN) for every particle in all the snap- 
shots. The defect fraction D (defined as the ratio of 
particles with a NN number different from 6 and the to- 
tal particle number) was then determined for every snap- 
shot. Using the neighbor information tables, the complex 
bond-order parameter Ge{j) — ^ ^xp [i69j(^)] was 
calculated for each particle, where j labels the particles, 
Nj is the number of true nearest neighbors of particle j, 
and Qj{l) is the angle of the l^^ nearest neighbor bond 
from an arbitrarily chosen, but fixed direction (the x-axis 
in our case) . Based on the spatial distribution of Ge , the 
visualization and identification of the crystallite structure 
is possible, since the magnitude jGel is close to unity in- 
side a grain but small for particles situated within the do- 
main walls (Fig.[IJb)). The complex argument Arg(G6), 
on the other hand, measures the overall angular orien- 
tation of a given particle neighborhood, or entire grains 
(Fig. [He)). Both the g{r) pair correlation function and 
the geir) = {Gl{0)GQ{r)) bond-order correlation func- 
tion were calculated for every snapshot. Assuming ex- 
ponential decay for the envelopes of both g{r) — 1 and 
geif) at large r, the corresponding correlation lengths 
and ^6 were then obtained employing an upper envelope 
least-square fit for each snapshot (Fig. [T](d,e)). To avoid 
boundary effects, both the correlation lengths and de- 
fect fractions were evaluated for particles situated within 
the red circles shown in Fig. [TJa) (which has a diameter 
0.9L). Distances appear normalized to the 2D Wigner- 
Seitz radius a = 1/ y/im, with n being the surface density 
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FIG. 1. (color online) Illustration of the data evaluation pro- 
cess for a liquid (left column) and a solid (right column) using 
snapshot data from experiment. Row (a): particle positions 
as identified in the raw images. Rows (b) and (c): magnitude 
and complex phase angle, respectively, of G6{r) bond order 
parameter. Rows (d) and (e): \g{r) — 1| (pair) and g6{r) 
(bond-order) correlation functions of the sample snapshots 
involving particles situated within the red circles in (a). 

of the detected particles. 

The time-evolution of ^6{t)i ^g{t) ^nd 1/D{i) (the in- 
verse defect fraction) is depicted in Fig. [2] for a single. 
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FIG. 2. (color online) Time evolution of the correlation 
lengths ^6 (red squares), (green circles) and the inverse 
defect fraction 1/D (blue triangles) for a single, representa- 
tive freezing cycle. Lines represent linear fits on a log- log plot 
for the short and long time domains. 

representative quenching cycle. One can easily identify 
two apparently linear regimes on the log-log plot in the 
case of the correlation functions, showing different time- 
evolution on the short and long time-scales. Similar ob- 
servations have been made for colloidal suspensions as 
reported qualitatively in [3] • In our case, the short time 
period lasts about 10-15 typical oscillation cycles for a 
single particle located in its equilibrium lattice position 
with a frequency i^dust- The defect density shows, on 
the other hand, a single characteristic decay time scale. 
This can be understood realizing, that: (i) at very early 
times after the temperature drop the particles ballisti- 
cally move towards their local equilibrium positions, and 
although crystallites are not yet formed, the defect den- 
sity is high and its distribution is homogenous, (ii) After 
some time, small crystallites fill out the space leaving 
defects to be accumulated in the walls between those 
grains. In this case (assuming self-similar structure of 
the domain walls at different times) the relationship be- 
tween the domain size and the inverse defect fraction can 
be assumed to be linear, (iii) Gg, and thus (?6(^) are 
continuous quantities, their values can increase for every 
particle with increasing order. The number of neighbors, 
on the other hand, is a discrete variable, taking on only 
integer values, (iv) Finally, the defect fraction "sensing" 
the particles with number of neighbors unequal to six, 
is less sensitive to the increase of the initially very small 
bond-angular and positional order at early times. 

Fitting linear functions to the log-log representation of 
the time-evolution data, as shown in Fig.[2l results in the 
a exponents. Averaging over 22 individual freezing cycles 
(of the same particle cloud) results in: ag (early) = 1.0± 
0.16, ae{\a.te) = 0.39 ±0.1 for the bond-order correlation 
length; ag(early) = 0.87±0.25, ag(late) = 0.34±0.18 for 
the translational correlation length; and an = 0.40±0.05 



for the inverse defect fraction. 

Although the a exponents could be derived using the 
arbitrary units related to the experimental conditions, 
in order to provide a better understanding of the sys- 
tem, we have derived the following physical quantities 
based on the dynamical single-particle and collective be- 
havior in the solid phase and reference wave dispersion 
data from Yukawa lattice calculations (i) Wigner- 

Seitz radius a = 400 ± 15 /im, (ii) Debye screening 
length Ad/o = 1-5 ± 0.3, (iii) nominal plasma frequency 
ujo = 55±5 rad/sec, (iv) particle charge Q/e = 10''±10%. 

Detailed investigation into the time evolution of the 
quenching process shows that the dust temperature 
(measured in terms of kinetic energy) drops from about 
200Ti„ to 0.4T„i, where Tm is the critical value associated 
with the equilibrium solid-liquid transition temperature 
[2^. This temperature drop can be well fitted with an 
exponential decay in the form T{t) oc exp(— i/r) with 
T = O.lliO.Ol s, which is very close to u'^ = 0.15±0.05 s 
the estimated inverse dust-neutral collision frequency. A 
rapid increase of the observed correlation lengths is ini- 
tiated immediately when T reaches Tm = f"'(i^^)m/2fcB 
in its descent, where {v'^)ni = WQa^/Fm, with Fm w 157 
for Yukawa interaction with Ad/q = 1-5 [1^. The dura- 
tion of this early period is equal to the time needed for a 
particle (moving with the observed average velocity) to 
travel a distance of w 2a, which is very close to the lattice 
constant in the final ground state hexagonal configura- 
tion, supporting the above described qualitative picture 
of rapid early ordering induced by single-particle motion. 

Our molecular dynamics simulations were carried out 
for a large 2D particle ensemble (with N = 1 154 000 
particles) , interacting through a Yukawa type interaction 
potential $ ~ exp(— t/Ad)/?', where the Debye screening 
length was chosen to be half the ground state hexago- 
nal lattice constant (Ad — 0.95a). The topical review 
article [28| provides an in-depth description of the com- 
putational methodology followed in the present work. 
Quenching in the simulation was achieved by suddenly 
setting particle velocities to zero and applying a strong 
frictional force for a short period of time. This allowed 
changes in the temperature from twice to 0.34 times the 
melting temperature in less then lOa;^^^ plasma oscilla- 
tion cycles. Processing of the simulation data was per- 
formed in the same manner as previously introduced for 
the experiment. These results are shown in Fig. [3] 

To identify the grains, a direct size measuring algo- 
rithm similar to the flood fill method was adopted. The 
technique is very common in graphical applications where 
it is often used to fill closed areas with a selected color. In 
our case we start "coloring" assuming 1600 initial points 
evenly distributed over our simulation cell. A particle is 
considered to belong to the same domain as the starting 
point if (i) it belongs to the neighborhood according to 
the Delaunay triangulation, and (ii) its complex phase of 
Gq is within the ±15 degrees range with respect to that 
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FIG. 3. (color online) Time evolution of the correlation 
lengths ^6, (,g, the inverse defect fraction 1/D, and the size 
measured with our "flood fill" algorithm from MD simula- 
tion. Time is normalized to the 2D nominal plasma frequency 
(jjQ — \J nQ^ l^e^ma. Points are shifted for clarity. 



often significant and irreproducible influence of the ac- 
tual conditions in a particular complex plasma experi- 
ment llj. This result is supported by our molecular 
dynamics simulation for a 2D Yukawa system resulting 
in (J^^ — 0.375 ± 0.03. Comparison of our findings re- 
veals satisfactory agreement with earlier results concern- 
ing systems with low friction (0.3 < a < 0.45) and with 
liquid (overdamped) environment {a < 0.3) exhibiting a 
somewhat faster coarsening 
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of the starting point. This algorithm results in the num- 
ber of particles (the weight) of the domain containing 
the starting point. In the calculation of the average do- 
main weight (VF) we neglect cases with < 6, as these 
represent cases where the starting point was accidentally 
picked in a domain wall. The average linear size of the 
domains is approximated as Rw oc \J (VF). 

Performing linear fitting to the data as shown in Fig. [3] 
for times tuj^ > 500 results in the following exponents: 
ae ~ 0.41, ttg « 0.35, an ~ 0.36, and aw ~ 0.38. The 
values obtained by these four independent measures pro- 
vide consistent results if we assume for the real time evo- 
lution of the linear domain size a power law time depen- 
dence with an exponent a = 0.375 ± 0.03, corresponding 
to (and extending) the long time characteristics of the 
experiment. 

In conclusion, we have performed a series of com- 
plex plasma experiments measuring the time evolution 
of the crystallization process of a single layer system 
after a rapid temperature quench from the liquid to 
solid state. Our experiments show a rapid early pro- 
cess, driven by single particle motion, where the system 
undergoes a continuous transition from its initial amor- 
phous state to a polycrystalline state of small domains. 
After about 10 - 15 particle oscillation cycles, pattern 
formation slows and further grain coarsening is only pos- 
sible by collective rearrangement, resulting in the merg- 
ing of small domains to form larger ones. The linear 
domain size shows a power-law type time evolution with 
an exponent a*^'^''^^ = 0.93 ± 0.1 at the short times, while 
Q,iate _ Q 3g _|_ g Qg \ong times. The long-time char- 
acteristics could be observed under stable temperature 
conditions, which made the real time-dependence study 
possible, as it is not altered by systematic changes in the 
temperature. It reveals the true nature of grain coars- 
ening and pattern formation in 2D, not affected by the 
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